Characterization of the Solutions Set of Least-squares Problems by an Extension of Kaczmarz's Projections Method

We give a new characterization of the solutions set of the general (inconsistent) linear least-squares problem using the set of limit-points of an extended version of the classical Kaczmarz's projections method.We also obtain a \ step error reduction formula" which,at least from a theoretical view point,can give us information about the convergence properties of the algorithm.In section 3 we present a modi ed version of the above extended algorithm,obtained by a class of nonorthogonal transformations of the initial least-squares problem and we prove for it similar results as for the initial one. Some comparisons concerning computational aspects,between our algorithms and others existent in the literature ,are made in the last section of the paper. 2 1.Introduction Let A be an m n real matrix (e.g. n m for simplifying the presentation) and b 2 IRm a given vector. Let At; R(A); N(A); r(A) be the transpose,range,null-space and rank of A respectively.If S IRq is a closed convex subset,by PS we shall denote the orthogonal projection on S (all these are considered with respect to the Euclidean scalar product and norm denoted by < ; > and k k respectively).If b 2 R(A) (1) we shall denote by S(A; b) the set of all (classical) solutions of the consistent system Ax = b: (2) In the general (inconsistent) case b = 2 IR(A) we shall consider the linear least-squares problem : nd x 2 IRn such that k Ax b k= minfk Ax b k; x 2 IRng (3) and we shall denote by LSS(A; b) the set of all its solutions.It is well known (see e.g.[3]) that LSS(A; b) contains a unique element (denoted here by xLS and called the minimal norm least-squares solution) such that k xLS k k x k; (8)x 2 LSS(A; b): (4) Let b2 2 IRm be de ned by b2 = PR(A)(b): (5) Then it can be easily proved that LSS(A; b) = S(A; b2): (6) There exist characterizations of the set LSS(A; b) based on the normal equation associated to (3),singular values or QR decompositions of A and generalized inverses of the matrix A (see e.g. [2],[3],[14]).In Section 2 of the present paper we shall present another possibility to characterize the elements of LSS(A; b),namely,as the set of all limit-points (with respect to the initial choice) of a sequence generated with an iterative method which extends in the general case of inconsistent problems (3) the classical Kaczmarz's projections method.(see [6],[10],[15]). This characterization is also 3 important because we obtain evaluations of the \step error reduction factor" of the algorithm and we give necessary and su cient conditions in order to obtain as limit-point the minimal norm solution xLS from (4). In Section 3 we prove that it is possible to introduce new projections in each iteration step of our algorithm ,without changing the set of its limit-points.This fact is important because (for some classes of problems) we can obtain in this way an acceleration of the convergence of the algorithm or,if a special set of new directions for projection are introduced we can obtain a direct solver for (general) inconsistent problems (3)-(4) (see e.g.[4],[8],[9],[11]).In Section 4 we brie y describe some other versions of iterative algorithms for approximating xLS in the general inconsistent case for (3) and we compare them with our algorithm from the view point of accuracy of approximation . In the rest of this introductory section we shall brie y present some basic facts about the (classical) Kaczmarz's projections method. We shall essentially follow in presentation the paper [15]. Let ai; i = 1; : : : ; m; j; j = 1; : : : ; n be the i-th row and the j-th column of the matrix A,respectively and bi the i-th component of the vector b. We shall suppose that ai 6= 0; j 6= 0; (8)i = 1; : : : ; m; j = 1; : : : ; n (7) and we shall consider the applications fi(b; ); F (b; ) : IRn ! IRn and 'j ; : IRm ! IRm de ned by fi(b; x) = x < x; ai > bi k ai k2 ai (8) F (b; x) = (f1 : : : fm)(b; x); x 2 IRn; (9) 'j(y) = y < y; j > k j k2 j ; j = 1; : : : ; n; (10) (y) = ('1 : : : 'n)(y); y 2 IRm: (11) Let also Pi; Q : IRn ! IRn be the linear applications (matrices) de ned by Pi(x) = x < x; ai > k ai k2 ai; (12) Q = P1P2 : : :Pm (13) 4 and R the real n m matrix of which the i-th column (denoted by (R)i) is (R)i = 1 k ai k2P1P2 : : :Pi 1(ai) (14) with P0= identity. Lemma 1.1(i) We have the equalities F (b; x) = Qx+ Rb;Q+ RA = I ; (15) (ii)N(A) and R(At) are invariant subspaces for Q;we have Q = PN(A) ~ Q (16) where ~ Q is the n n matrix de ned by ~ Q = Q PR(At); (17) (iii)The matrix ~ Q has the property k ~ Q k= supfk Qx k; x 2 R(At); k x k= 1g < 1; (18) (iv) N(At) and R(A) are invariant subspaces for the application ;we have = PN(At) ~ ; (19) where ~ is de ned by ~ = PR(A); (20) (v)The matrix ~ has the property k ~ k= supfk (y) k; y 2 R(A); k y k= 1g < 1; (21) (vi)For any y 2 IRm the matrix R from (14) satis es Ry 2 R(At): (22) Lemma 1.2.Let (xk)k 0 IRn; (yk)k 0 IRm be the sequences de ned by xk+1 = F (b; xk); k 0; x0 2 IRn; (23) yk+1 = (yk); k 0; y0 = b 2 IRm; (24) 5 and G the n m matrix G = (I ~ Q) 1R: (25) Then lim k!1 xk = PN(A)(x0) + Gb; (26) lim k!1 yk = PN(At)(b) = b1: (27) The matrix G satis es (GA)t = GA (28) and b 2 R(A)() AGb = b: (29) Remark 1.1. For the proofs of the above results see [15] and [10].. Remark 1.2. The iteration (23) is the classical Kaczmarz's projections method.(see [15]). 2.Characterization of the set LSS(A;b) We shall denote by LK(A; b) the set of all limit-points of the (classical) Kaczmarz's sequence (23) (with respect to the initial choice x0 2 IRn),i.e. (see (26)) LK(A; b) = fPN(A)(x0) +Gb; x0 2 IRng: (30) The following result gives,in the particular (consistent) case b = b2 2 R(A); (31) a rst characterization of the set LSS(A; b) using the limit-points of Kaczmarz's sequence (23). Proposition 2.1.For any matrix A satisfying (7),if the vector b 2 IRm is as in (31) then LSS(A; b) = LK(A; b) (32) and the minimal norm solution xLS is obtained as a limit point in (23) if and only if x0 2 R(At).In this case we have xLS = Gb (33) 6 Proof.If x 2 LK(A; b) then x = PN(A)(x0) + Gb for some x0 2 IRn and,using (31),(6) and (29) we obtain Ax = A(PN(A)(x0) + Gb) = AGb = b,thus x 2 LSS(A; b).Reversely,if x is an element of LSS(A; b),using again (6) and (31) we get Ax = b.If we de ne the vector x0 2 IRn by x0 = x Gb (34) we obtain,using (29),Ax0 = Ax AGb = b b = 0 i.e. x0 2 N(A), thus x0 = PN(A)(x0).From this last equality,(34) and the de nition (30) it holds that x 2 LSS(A; b) and thus equality (32). For the second part of the proof let x 2 LSS(A; b) and z 2 N(A) be two arbitrary vectors.Using (6) and (28) we obtain < z;Gb >=< z;GAx >=< Az;Gtx >= 0: (35) But,from (30) and (32),it holds that x = PN(A)(x0) + Gb, for some x0 2 IRn.Thus,using (35) with z = PN(A)(x0) and the orthogonal decomposition IRn = N(A) R(At) we obtain k x k2=k PN(A)(x0) k2 + k Gb k2 k Gb k2 : (36) The inequality (36) becomes equality if and only if PN(A)(x0) = 0 which is equivalent with x0 2 R(At).Then,the equality (33) obviously holds and the proof is complete. As we have already mentioned,the above result o ers a characterization of the set LSS(A; b) only in the particular case (31).For the general case of inconsistent problems (3) we shall use the following algorithm ,called (KE) (Kaczmarz Extended),which we have proposed in a particular form in [10]. Algorithm (KE) Let y0 = b; x0 2 IRn (arbitrary) and xk 2 IRn an already computed approximation.The next one xk+1 is obtained by yk+1 = (yk); (37) bk+1 = b yk+1; (38) xk+1 = F (bk+1; xk); k 0: (39) Remark 2.1. In the paper [10] we proved the convergence of the above de ned sequence (xk)k 0 to xLS only in the particular case x0 = 0.In what 7 follows,by a deeper analysis,we shall obtain the set of all limit-points of the algorithm (KE) and we shall use it for the characterization of LSS(A; b) in the general inconsistent case.The following result gives a rst information about the sequence (xk)k 0 for an arbitrary choice x0 2 IRn. Lemma 2.2.Let (xk)k 0 be the sequence generated with the algorithm (KE),with an arbitrary initial choice x0 2 IRn.Then PN(A)(xk) = PN(A)(x0); (8)k 0: (40) Proof.Let k 0 be arbitrary xed and suppose that (for this k) we have PN(A)(xk) = PN(A)(x0): (41) >From (15),(16),(37)-(39) and (41) we obtain xk+1 = Qxk + Rbk+1 = PN(A)(x0) + ~ Qxk +Rbk+1: (42) Using (17) and the fact that R(At) is an invariant subspace for Q we obtain ~ Qxk 2 R(At) (43) which together with Rbk+1 2 R(At) (see (22)) gives us ~ Qxk + Rbk+1 2 R(At): (44) >From (42),(44) and the decomposition IRn = N(A) R(At) we obtain that PN(A)(xk+1) = PN(A)(x0); (45) and the proof is complete. Proposition 2.3.For any matrix A satisfying (7) and any vector b 2 IRm the sequence generated with the algorithm (KE) is convergent and lim k !1 xk = PN(A)(x0) +Gb2: (46) Proof.Using (37)-(39),(5),(27),(15),(25) and (42), we succesively obtain xk+1 (PN(A)(x0) +Gb2) = ~ Qxk +Rbk+1 Gb2 = = ~ Qxk + Rbk+1 [(I ~ Q) + ~ Q](I ~ Q) 1Rb2 = (47) = ( ~ Qxk ~ QGb2)+(Rbk+1 Rb2) = ~ Q[xk (PN(A)(x0)+Gb2)]+R(b1 yk+1): 8 Then,as in [10],using (27),(24),(37) and (19) we get b1 yk+1 = b1 k+1b = b1 (PN(At) ~ k+1)b = ~ k+1b; (48) thus,following (47) xk+1 (PN(A)(x0) +Gb2) = ~ Q[xk (PN(A)(x0) +Gb2)] R~ k+1b; (8)k 0: (49) Now,if we denote by ek the error vector ek = xk (PN(A)(x0) + Gb2); k 0; (50) a recursive application of (49) gives us ek = ~ Qke0 k 1 X j=1 ~ Qk jR~ jb R~ kb; k 2; (51) thus,by tacking norms, k ek k k ~ Q kkk e0 k + k R kk b k (ck+ k ~ kk); k 2; (52) where by (ck)k 2 we denoted the sequence of positive numbers de ned by ck = k 1 X j=1 k ~ Q kk jk ~ kj ; k 2: (53) Let (ak)k 1; (bk)k 1 be the sequences of positive numbers de ned by ak =k ~ kk; bk =k ~ Q kk; k 1: (54) By observing that (see also (18) and (21)) lim k !1 ak = 0; (55) 1 X k=1 bk <1; (56) and using a similar argument as in [10](theorem 3.1) we obtain lim k !1 ck = 0; (57) 9 which together with (52) (and using again (18) and (21)) gives us lim k !1 k ek k= 0 (58) and the proof is complete. Remark 2.2.From (49) and (50) we obtain the inequality k ek+1 k k ~ Q kk ek k + k R kk b k k ~ kk+1; k 0; (59) which tells us that the convergence properties of the algorithm (KE) are \controled" by the values k ~ Q k and k ~ k (from (18) and (21)).This fact can be a starting point in the analysis of the possibility of accelerating the convergence properties of the algorithm (KE) (see in this sense [4],[9]).More than that,in the particular case (31),from (27) we obtain that the sequence (xk)k 0 generated with the algorithm (KE) is the classical Kaczmarz's sequence (23) (more clear,we have yk = 0; (8)k 1).Then,from (49) we obtain the following relations xk+1 (PN(A)(x0) +Gb) = ~ Q[xk (PN(A)(x0) +Gb)] (60) and k ek+1 k k ~ Q kk ek k; k 0; (61) which tells us that,if (31) holds,then the convergence properties of the classical Kaczmarz's method (23) are \controled" by k ~ Q k (in fact,in this case k ~ Q k becomes a \step error reduction factor",see e.g. [16]). Remark 2.3.If we denote by LKE(A; b) the set of all limit-points of the algorithm (KE) (with respect to the initial choice x0 2 IRn) from Proposition 2.3 we get LKE(A; b) = fPN(A)(x0) +Gb2; x0 2 IRng: (62) Now we can prove our main result,announced at the begining of this section (similar with Proposition 2.1,but in the general inconsistent case). Theorem 2.4.(i) For any matrix A satisfying (7) and any vector b 2 IRm we have the equality LKE(A; b) = LSS(A; b) (63) (ii) The minimal norm solution xLS is obtained as a limit-point in (37)-(39) if and only if x0 2 R(At). Proof.From (6),(30) and Proposition 2.1 we get the equality LSS(A; b) = LSS(A; b2) = fPN(A)(x0) + Gb2; x0 2 IRng: (64) 10 This equality together with (62) completes our proof (the assertion (ii) is obvious). Remark 2.4.The condition x0 2 R(At) can be easily ful led.Indeed, it su ces to take x0 = 0 or x0 = ai (i.e. a row of A). 3.Invariance of the set LKE(A;b) under a class of non-orthogonal transformations It is well known that the set LSS(A; b) is invariant under orthogonal transformations of the problem (3).In fact,this is the basic idea used in solving (3) by direct methods (as complete QR decomposition,see e.g. [3]).It is also very clear that,in general, non-orthogonal transformations change the set LSS(A; b) and also the minimal norm solution xLS (even in the case of transformations using diagonal matrices;see e.g.[5]). In this section we shall prove that it exists a class of non-orthogonal linear transformations of (3) for which the set of limit-points of a corresponding \transformed" algorithm of the type (KE), is too LKE(A; b).The class of the above mentioned non-orthogonal transformations is considered following the idea of \supplementary directions for projection" used by the author in [9] (see also [4]),in steps (37) and (39),in order to obtain an improvement of the convergence properties of the algorithm (KE) or a direct solver (see [11]).For simplicity we shall present the construction of the \transformed" algorithm and the corresponding convergence results only for one new direction in each of the steps (37) and (39).After the presentation it will be clear that a simple recursive argument extend the construction and results to an arbitrary ( nite) number of new directions. Thus let a0 2 IRn; 0 2 IRm be two vectors obtained as linear combinations of the rows and columns of A,respectively,i.e. a0 = m Xi=1 iai; 0 = n Xj=1 j j ; (65) with i; j 2 IR.Let also b0 2 IR be the new component of b de ned by (see (65)) b0 = m Xi=1 ibi: (66) 11 We suppose that (see also (7))a0 6= 0; 0 6= 0 (67) and let i 2 f1; : : : ; mg; j 2 f1; : : : ; ng be two arbitrary xed indices.We shall de ne the new projections f0(b; ) : IRn ! IRn and '0 : IRm ! IRm by (see (8),(10)) f0(b; x) = x < x; a0 > b0 k a0 k2 a0; (68) '0(y) = y < y; 0 > k 0 k2 0; (69) and the new applications F0(b; ) : IRn ! IRn and 0 : IRm ! IRm (see (9),(11)) F0(b; x) = (f1 : : : fi f0 fi+1 : : : fm)(b; x); (70) 0(y) = ('1 : : : 'j '0 'j+1 : : : 'n)(y): (71) Remark 3.1. For i = m we understand that in (70) projection f0 is applied rst,before fm (and similar for j = n in (71)). The \transformed" algorithm (of the type (37)-(39)) denoted by KE0 is the following. Algorithm (KE0) Let y0 = b; x0 2 IRn and xk 2 IRn an already computed approximation.The next one,xk+1 2 IRn,is obtained by yk+1 = 0(yk); (72) bk+1 = b yk+1; (73) xk+1 = F0(bk+1; xk): (74) Now,we shall write the above steps (72) and (74) in terms of some transformations of the matrix A and the vector b.For this,let rstly A be the (m+ 1) n matrix with the rows a1; : : : ; ai; a0; ai+1; : : : ; am (in this order) and b 2 IRm+1 the vector b = (b1; : : : ; bi; b0; bi+1; : : : ; bm): (75) 12 Let F ( b; ) : IRn ! IRn be the application de ned as in (9) but with A; b instead of A and b,i.e. F ( b; x) = ( f1 : : : fm+1)( b; x); (76) where fk( b; ) : IRn ! IRn are given by fk( b; x) = fk(b; x); k = 1; 2; : : : ; i fi+1( b; x) = f0(b; x) (77) fk( b; x) = fk 1(b; x); k = i+ 2; : : : ; m+ 1: Now,let  be them (n+1)matrix with columns 1; : : : ; i; 0; i+1; : : : ; m (in this order) and ̂ : IRm ! IRm the application de ned as in (11),but with  instead of A,i.e. ̂(y) = ('̂1 : : : '̂n+1)(y); (78) where '̂k : IRm ! IRm are given by '̂k(y) = 'k(y); k = 1; 2; : : : ; j '̂j+1(y) = '0(y); (79) '̂k(y) = 'k 1(y); k = j + 2; : : : ; n+ 1: It is very clear that we have the equalities F ( b; x) = F0(b; x); (8)x 2 IRn; (80) ̂(y) = 0(y); (8)y 2 IRm: (81) Then,we can write the steps (72) and (74) as follows yk+1 = ̂(yk); (82) xk+1 = F ( bk+1; xk): (83) Remark 3.2.By bk+1 in (83) we understand the vector from IRm+1 obtained from bk following (66) and (75). Remark 3.3.Let T be the (m+ 1) m matrix with the rows 1; : : : ; i; ; i+1; : : : ; m,in this order,where 1; : : : ; m are the vectors from the cannonical basis of IRm,i.e. k = (0; : : : ; 1; : : : ; 0); (84) 13 with 1 on position k and 2 IRm is given by (see (65)-(66)) = ( 1; 2; : : : ; m): (85) Then we obviously have A = TA; b = Tb: (86) Remark 3.4.We replay the above construction using At instead of A (i.e. 1; 2; : : : ; n are the rows of At) and we obtain an (n + 1) n matrix denoted by T1 with the propertyÂt = T1At: (87) Then,by de ning the n (n+ 1) matrix by = T t 1 (88) we obtain the equality  = A : (89) The above remarks tell us that steps (72) and (74) of the new algorithm (KE0) are constructed accordingly to the linear (nonorthogonal) transformations (86) and (89) of the matrix A and the vector b. Now let Q; ~ Q; R; G be the matrices constructed as in (13),(17),(14),(25),respectively,but with A instead of A, and ̂; ~̂ constructed as in (13),(20),respectively but with  instead of A.Then,all the results of sections 1 and 2 rest true for these applications (with respect to the corresponding elements ,i.e. A; b and Â).More than that we have Lemma 3.1. The following are true : N(Ât) = N(At); (90) Tyk+1 = T (b1 + ~̂ k+1b); (91) with yk+1 from (73). Proof.Using the fact that 0 is a linear combination of 1; : : : ; n (see (65)) we have the following sequence of equivalences proving (90) Âty = 0, f< j ; y >= 0; j = 1; : : : ; n; < 0; y >= 0g , , f< j ; y >= 0; j = 1; : : : ; ng , Aty = 0: 14 >From (19) and (24) (with  instead of A) we obtain Tyk+1 = T (̂k+1b) = T (PN(Ât)(b) + ~̂ k+1b): (92) But,from (90) and (27) we obtain PN(Ât)(b) = PN(At)(b) = b1 (93) which together with (92) gives (91) and the proof is complete. Now we can prove the following result concerning our algorithm(KE0). Proposition 3.2.The sequence (xk)k 0 generated with the algorithm (KE0) is convergent and we have lim k !1 xk = PN( A)(x0) + GTb2: (94) Proof.As in the proof of Proposition 2.3 we obtain (see (47)) xk+1 (PN( A)(x0)+ GTb2) = ~ Q[xk (PN( A)(x0)+ GTb2)]+ RTbk+1 RTb2: (95) Using (91) and (93) we also obtain (see (73)) RTbk+1 RTb2 = RTb1 RTyk+1 = RT ~̂ k+1b which together with (95) gives us xk+1 (PN( A)(x0)+ GTb2) = ~ Q[xk (PN( A)(x0)+ GTb2)] RT ~̂ k+1b; (8)k 0: (96) Then,because from (18) and (21) we have k ~ Q k< 1; k ~̂ k< 1; (97) by using similar arguments as in (50)-(58) we can complete the proof. Let now LKE0(A; b) be the set of all limit-points of the algorithm (KE0) (with respect to the initial approximation x0).The main result of this section is the following. Theorem 3.3.We have the equality LKE0(A; b) = LKE(A; b): (98) Proof.Using a similar argument as in the proof of (90) we obtain N( A) = N(A): (99) 15 Then,following the de nition of A we obtain S(A; b2) = S( A;Tb2) (100) from which,using similar arguments as in the proof of Proposition 2.1 we get Gb2 = GTb2: (101) >From (89),(91),(94) and the de nition of the set LKE0(A; b) we obtain LKE0(A; b) = fPN(A)(x0) +Gb2; x0 2 IRng; (102) which together with (62)-(63) gives us (98) and the proof is complete. 4.Remarks concerning computational aspects Unfortunately the convergence properties of the algorithm (KE) (or its modi ed version (KE0)0 are ,generally, bad.This because they use in steps (37) and (39) (respectively, (72) and (74)) succesive projections similar with the classical Kaczmarz's iteration (23) (for which is well known that has slow convergence).But,beside this bad thing,it exists a \good" computational property of the above mentioned algorithms : their \compact" form,i.e. a single iterative step which is repeated until a desired accuracy (with respect only to the sequence (xk)k 0) is achieved.We can use di erent \stopping tests" for these algorithms which can be easily programmed and which ensure only one control at the end of each iteration,as e.g. one of the following (" is the desired accuracy) : k xk+1 xk k "; (103) maxfjxk+1 i xki j; i = 1; : : : ; ng "; (104) or (see e.g. [10]) k Axk bk k ": (105) (with bk from (38)). Concerning other iterative algorithms for general inconsistent problems as (3),at least at author's knowledge,there exist in the corresponding literature the following ones A. An algorithm due to L.D.Pyle (see [13]) which is composed by the following three di erent iterative steps : Step 1. Compute b2 = lim k!1(b yk) (106) 16 where the sequence (yk)k 0 is generated as in (37) with y0 = b. Step 2. Compute an arbitrary solution x̂ 2 IRn of the system Ax̂ = b2 (107) This is done by the author in (1 q + 1) applications of an iterative procedure similar with the classical Kaczmarz (23).The number q 0 is the dimension of N(A) and the value of can not be apriori determined (it depends on the positions of lineary independent columns ofA in the matrix). Step 3. Compute xLS = lim k!1(x̂ xk); (108) where the sequence (xk)k 0 is generated by x0 = x̂; xk+1 = f(xk); k 0; (109) where f : IRn ! IRn is an application de ned as from (10)-(11),but using the rows ai of A instead of its columns j . B. An algorithm due to R.J.Plemmons (see [1]) which is composed by two di erent iterative steps as follows: Step 1.Compute an arbitrary least-squares solution x of (3) using a conjugate gradient-like method (CGPCNE in [1]) applied to the normal equation of (3),i.e. AtAx = Atb: (110) Then,compute the corresponding residual r and the vector b2 (from (5)) by r = b Ax ; b2 = b r : (111) Step 2. Compute the minimal norm solution of the (consistent) system Ax = b2 (112) (which is the desired xLS ;see also (6)),by another type of conjugate gradient algorithm (CGPCMN in [1]). C. An algorithm proposed (and used) by the author in [8],composed by the following two separated iterative steps : Step 1. Compute the vector b1 from (27) using (37) with y0 = b. Step 2. Compute xLS by applying the classical Kaczmarz's iteration (23) (with x0 = 0) to the (consistent) system Ax = b b1: (113) 17 All the above algorithms have the same \big" de ect : in Step 2 only an ap-proximation of the requested value (or vector) comes from Step 1,a.s.o.Thus,atthe end , a superpositioning of approximations can appear,which (at leastfor ill-conditioned matrices A) can generate serious troubles.We shall ilus-trate this in what follows by a simple example involving the above algorithmC,but similar computations can be performed also with respect to A andB.Let A = (aij)i;j be the (n+ 1) n real matrix given byaii = 2; i = 1; : : : ; n 1; ai;i 1 = ai;i+1 = 1; i = 2; : : : ; n 2; (114)a12 = an 1;n 2 = an2 = an+1;1 = 1; a1n = an 1;n = an1 = an+1;2 = 1;and b = (bi)i 2 IRn+1 the vectorb1 = bn 1 = 2; bn+1 = 10(115)(the other entries of A and b are equal to zero).From construction r(A) =n 1; b =2 R(A) i.e. we are in the general case of an inconsistent and rank-de cient problem (3) (see e.g. [3]).We applied the above algorithm C asfollows :starting with y0 = b,compute with (37) a vector y such thatk At y k "y ;(116)starting with x0 = 0 and using (23) with b = b y compute the rstapproximation xk such thatk xk+1 xk k 10 6:(117)Then we tried to obtain the approximation xk with our algorithm (KE),i.e.starting with x0 = 0; y0 = b we tried to compute an xq such thatk xq xk k 10 6:(118)We rstly observe that all the above computations depend on the accuracy"y with which we approximated b1 in the rst step of C (see (116).In table 1we present the bigest values of "y for which the approximation in (118) waspossible (for di erent values of n).This means that, for example,in the casen = 15,for "y = 10 5 we could not get the corresponding xk with a term ofthe sequence(xq)q 0.But we already know that (see Theorem 2.4)limq!1 xq = xLS :(119)18 Thus,we conclude that in this case the value xk is \too far" from xLS (atleast at a distance > 10 6).In order to overcome this di culty we have tochoose a smaller value for "y in (116),and, unfortunately,we can observe intable 1 that this depends (at least) on the dimensions of the matrix A (i.e.bigger n is , smaller must be the accuracy "y). Thus,it is possible that forvery big values of n we need too small values for "y which will deterioratethe computations.Table 1.n maximal value of "y610 61510 72510 83510 8References[1] Bjork A.,Elfving T. Accelerated projection methods for computingpseudoinverse solutions of systems of linear equations; BIT 19(1979),pp.145-163.[2] Boullion T.L.,Odell P.L. Generalized inverse matrices;WilleyInter-science,New York,1971.[3] Golub G.H.,van Loan C.F. Matrix computations;The John's Hop-kins Univ.Press,Baltimore,1983.[4] Golubovici G.,Popa C. Supplementary relaxations for the accelera-tion of a class of iterative methods;Analele Stiinti ce ale Univ."OVIDIUS"Constanta,Seria Matematica,vol.III,fasc.1 (1995),pp.52-62.[5] Hanke M.,NiethammerW. On the acceleration of Kaczmarz's methodfor inconsistent linear systems;Linear Alg.Appl., 130 (1990),pp.83-98.[6] Kaczmarz S. Angenaherte Au osung von Systemen linearer Gle-ichungen;Bull.Acad.Polonaise Sci. et Lettres A (1937), pp.355-357.[7] Mc Cormick S.,Ruge J. Unigrid for multigrid simulation; Math. ofComp.,41 (1983),pp.43-62.19 [8] Popa C. An iterative method for CVBEM systems.Part I : TheKaczmarz algorithm;Part II:The unigrid method;Advances in EngineeringSoftware,16 (1993),pp.61-69.[9] Popa C. Directional relaxations for the acceleration of Kaczmarz'siterative method;Preprint,Rechenzentrum Univ. Karlsruhe,Germany,June1994.[10] Popa C. Least-Squares Solution of Overdetermined InconsistentLinear Systems Using Kaczmarz's Relaxation;Intern.J.Computer Math., 55(1995),pp.79-89.[11] Popa C. Succesive orthogonal projections for exact computationof minimal norm solution of inconsistent and rank-de cient least-squaresproblems;Preprint,CS96-19, December 1996,Dept.Appl.Math.Comp.Sci.,WeizmannInstitute of Sci.,Israel.[12] Popa C. An extension of Kaczmarz's projections method with relax-ation parameter to inconsistent and rank-de cient problems;Preprint,CS96-14,October 1996,Dept.Appl.Math.Comp.Sci.,Weizmann Institute of Sci.,Israel;submitted to BIT.[13] Pyle L.D. A generalised inverse -algorithm for constructing inter-section projection matrices with applications;Numer.Math.,10 (1967),pp.86-102.[14] Radhakrishna R.,Mitra S.K. Generalized inverse of matrices andits applications;John Willey and sons Inc.,New York,1971.[15] Tanabe K. Projection Method for Solving a Singular System ofLinear Equations and its Applications;Numer.Math., 17(1971),pp.203-214.[16] Young D.M. Iterative solution of large linear systems. AcademicPress,New York,1971.20